Global theta:Int=45


Global gamma2:Double=0.3*0.3
Global sigma2:Double=2.8*2.8 '4.5*4.5 '2.8*2.8

Global lambda:Double=3.5 '5.6 '3.5

Global w:Int=3



Function gabor:Double(x:Int, y:Int, theta_:Int)
  Local xx:Double=Double x*Cos(theta_) + Double y*Sin(theta_)
  Local yy:Double=Double -x*Sin(theta_) + Double y*Cos(theta_)

  Return Exp(-(xx*xx + gamma2*yy*yy) / (2*sigma2)) * Cos(Double xx * Double 360/lambda)
End Function






Graphics 800, 600


Global str:String="__constant__ float gabor"+theta+"["+(2*w+1)+"]["+(2*w+1)+"] = { "


Global tmp:Double=0
Global tmp2:Double=0



Local p:Double=0
For y=-w To w
  For x=-w To w
    Local f2:Float=gabor(x, y, theta)
    p=p+f2*f2
  Next
Next

p=Sqr(p)



Global sum:Double=0

For y=-w To w
  For x=-w To w
    Local x_:String=x
    Local y_:String=y
    If x>=0 Then x_="+"+x
    If y>=0 Then y_="+"+y

    Local f:Float=gabor(x, y, theta)
    f=f/p


    sum=sum+f
  Next
Next

sum=-sum/((2*w+1) * (2*w+1))






For y=-w To w
  str=str+"{"

  For x=-w To w
    x_=x
    y_=y
    If x>=0 Then x_="+"+x
    If y>=0 Then y_="+"+y

    f=getgabor( gabor(x, y, theta), p, sumpos, sumneg)


    tmp=tmp+f
    tmp2=tmp2+f*f


If mask(x, y, theta)
  str=str+f+", "
EndIf

    f=f+1
    SetColor f*128, f*128, f*128
    For i=0 To 2
      For j=0 To 2
         Plot (x+w)*2+i, (y+w)*2+j
      Next
    Next
    'Plot x+w, y+w
  Next

  str=str[0..str.length-2]+"}, "
Next


str=str[0..str.length-2]+" };"


Print tmp
Print tmp2
Print "-----"
Print


Print str




Flip
WaitKey
End




Function getgabor:Double(f:Double, p:Double, pos:Double, neg:Double)
  f=f/p+sum

Rem
  If f>0 Then f=f/pos
  If f<0 Then f=f/neg
End Rem

  Return f
End Function



Function mask(x, y, theta)
  Return 1

Rem
  If theta=45 Or theta=135
    If (x=-w And y=-w) Or (x=w And y=w) Or (x=-w And y=w) Or (x=w And y=-w)
      Return 0
    Else
      Return 1
    EndIf
  Else
    Return 1
  EndIf
End Rem
End Function





